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. In this paper we present an approach in which synchronization of chaos is used to address iden- 

O ' tification problems. In particular, we are able to identify: (i) the discontinuity points of systems 

■ described by piecewise dynamical equations and (ii) the delays of systems described by delay dif- 

ferential equations. Delays and discontinuities are widespread features of the dynamics of both 
natural and manmade systems. The foremost goal of the paper is to present a general and flexible 

43' 

Oh methodology that can be used in a broad variety of identification problems. 
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Recent work [1|, |2j, |3| has shown that synchronization of chaos can be conveniently used as a tool to identifying the 
dynamical equations of unknown systems |lfl| . For instance, in Q a largely unknown chaotic (nonlinear) system was 
, coupled to a model system and a general adaptive strategy was proposed to make them synchronize by adaptively 
varying the parameters of the model until they converge on those of the true system. The strategy takes advantage 

^" : n 

of the fact that complete synchronization of chaos is only possible if the coupled systems are exactly the same [2CJ ■ 
' In linear systems, an observer is a model dynamical system that is able to reconstruct the state of an unknown true 

o : n 

system from knowledge of its dynamical inputs and outputs. References [J] have outlined the connection between the 
problem of synchronization of dynamical systems and the problem of the design of an observer to reconstruct the 



state of an unknown system. The model system introduced in [2[ acts as an observer that dynamically identifies the 
parameters and the initial conditions of the nonlinear functions describing the dynamics of the true system, when 
they belong to a certain class (e.g., they are smooth and polynomial up to a given degree). In general, a condition 
that needs to be satisfied in order for the observer to fully reconstruct the state of the unknown system is that the 
inputs are persistently exciting signals, i.e, they solicit all the dynamical modes of the true system. In the case of 
chaotic autonomous systems, this requirement is naturally satisfied by the dynamics, which is in a state of persistent 
excitation, even in the absence of inputs [21] . This has important consequences since it allows to possibly extract a 
large number of unknown quantities from one chaotic time trace, without the need of introducing external inputs (see 



e.g. 2; la])- On the other hand, a limitation is that stability of synchronization of chaotic systems is usually evaluated 
locally about the synchronized evolution and for the parameters of the systems being the same (i.e., it relies on the 
successfulness of the identification strategy). Therefore, the effectiveness of the strategy may be sensitive to the choice 
of the initial conditions and it may depend on a careful selection of the adaptive strategy parameters. 

Chaos can also arise in systems whose dynamics is described by delay dynamical equations or by piecewise lin- 
ear/nonlinear equations. These are both very common situations in nature and in applications. For example, piecewise 



dynamical equations appear in mechanical systems with impacts 6J, relay-feedback systems 



jj, DC/DC converters 



8| . Delay dynamical equations are usually invoked to describe physiological processes in hematology, cardiology, neu 



rology, and psychiatry 9] ; but find application also in chemistry, engineering and technological systems 



are lasers subject to optical feedback 
14j . and traffic flow models 



12], mechanical vibrations 



ldl | ; examples 



131 ]. control engineering 



111 ], high-speed machining 
151 ] . Our goal is introducing a flexible adaptive strategy that can be used to identify 
either discontinuity points or delays of the equations of unknown chaotic systems from knowledge of their state evolu- 
tion. To our knowledge, these two important and interesting problems have not been addressed yet in the literature. 
More in general, the paper presents a general methodology that can be used to solve a broad variety of identification 
problems, including apparently difficult ones. 

As a first problem of interest, we consider that the unknown true system evolves in discrete time and is described 
by a piecewise dynamical equation. In particular, we are concerned with the case that the unknown system evolution 
obeys, 



x k+1 = U(x k ), 



(1) 



where x k £ R, f a (x) 6 R — > R is a piecewise-nonlinear function, 



U(x) = 



f^x), if g(x)<a, 



(2) 



otherwise, 



where f(x) E R -> R, f n (x) 6 R -> R, g{x) e R -> R and a is a scalar. For the sake of simplicity and without loss 
of generality, here we assume that the true system state a; is a scalar quantity. The more general case that the state 
is n-dimensional is discussed in an example that we present successively. 

Our goal is to identify the unknown parameter a from knowledge of the temporal evolution of x k . To this aim, we 



3 



try to model the dynamics of the true system by 

y k+1 = U(y k ), (3) 

where a' is an estimate of the unknown true coefficient a, and y is defined as, 

y k = ey k + (1 - e)x\ (4) 

e 6 R. Note that the model is coupled to the true system through x k in (|4]). We proceed under the assumption that 
the value of e in ((J) is such that when a' = er, the synchronized solution y — x, is stable. Our strategy (to be specified 
in what follows) seeks to synchronize the model with the true systems, by dynamically adjusting a' to match the 
unknown true value of a. We now introduce the following potential [16[, 

V = [x-y] 2 , (5) 

^ > by definition; ^ = when y — x, that is, when the true system and the model system are synchronized. Our 
goal is to minimize the potential (i.e., to achieve synchronization between the two systems) by dynamically adjusting 
the estimate a' . Thus, we propose to adaptively evolve a' in time, according to the following gradient descent relation, 

«' k+1 °' k = -v^p = M* k -y k fS = - y k ]v\ (6) 

?/ > 0. We note that the term p k = dy k /da' appears in Eq. ([B]); therefore, we seek to find a recurrence equation that 
describes the evolution of p k in time. We note that ([3]) can be rewritten as, 

yk+ l = f I { yk )H{(T > _ g{ y k)) + fl { f ){l U{p < _ g{ f ))]: (J) 
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where TL is the Heaviside step function. H(x) = 1 if x > 0, otherwise. This allows us to write, 



P 



k+l 



a k p k + b k , 



with, 



(8a) 



a k = e{Df{y k )H{a' - g(y k )) + Df H (y k )[l H{o> g(y k ))\ + [f n (y k ) - f(y k )]5{a' - g(y k ))Dg(y k )} : (8b) 



where S(-) is the Dirac delta function. ([8£i) is an auxiliary difference equation that completes the formulation of our 
adaptive strategy. In fact, by iterating the set of equations (|3|6|8|) and by using knowledge of the time-evolution of the 
true system state Xk, we obtain a time evolving estimate a 1 of the unknown quantity a. Thus, the adaptive strategy 
is fully described by the set of equations (|3I6I8|) . 

Figure 1 shows the results of a numerical experiment in which Eq. ^ is the tent map equation, f I (x) = fix, 
f II (x) = (1 — fi)x, fi — 1.4, a = 0.6. For these values of the parameters, the tent map dynamics is chaotic. Moreover, 
we choose g(x) = x, h(x) = x, e = 0.55. One of the definitions of the Dirac delta function is the following, 



In our numerical simulations in Fig. 1, we have approximated 5{u) by <5 e (u), with e = 0.1. We initialize both the 
true and the model systems from uniformly distributed random initial conditions between and 1, and a 1 is evolved 
from an initial estimate that is far away from the true value of tr, i.e., cr'(0) = 0.31. As can be seen from Fig. 1, 
our adaptive strategy is successful in identifying a and in achieving synchronization between the true and the model 
systems. 

Our approach can be extended to the case of continuous-time piecewise dynamical systems, and whose state is 
in general n-dimensional. Moreover, our strategy can be generalized to identify more than one discontinuity point. 
Though we do not report here a formal formulation of our strategy to encompass all these possible variations, in 
what follows we present a numerical example, which addresses the case of a continuous-time n-dimensional piecewise 
system, n — 3. Namely, we consider that the true system dynamics is described by the Chua equation, x(£) = 



& fc = [/V)-/"(y fc )]%' -git)) 



(8c) 



5{u) = lim 8 e (u) — 2 lim e 1 max(l — |w/e|,0) 




(9) 
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FIG. 1: (a) and (b) show the time evolutions of x(k) (in black, thin line) and y(k) (in grey, thick line) respectively at the 
beginning (k G [0, 50]) and at the end (k € [149950, 150000]) of the run (in (b) the two curves are superposed). Plot (c) shows 
\x(k) — y(k)\ versus k £ [0, 1.5 x 10 5 ]. As can be seen, perfect synchronization is achieved for k approximately greater than 
7.5 x 10 4 . (d) shows the time evolution of a'{k) (in grey), which converges to the true value of a = 0.6 (the a-th ordinate line 
is plotted as a black solid line in (d)). f I (x) = fix, f 11 (x) — /i(l — x), fi = 1.4, a = 0.6, g(x) — x, h(x) = x, e = 0.55, r\ — 10~ 5 , 

£ = 0.1. 



[x 1 (t),x 2 (t),x 3 (t)] T , x(t) = i^(x(t)), 



^(x) 



a(x 2 (t) - Xl (t) -MM*))> 
xi(t) - x 2 {t) + x 3 (t), 
-0x 2 (t), 



(10) 



where the piecewise scalar function 4> a {x) is defined as, 



0(7 (z) = < 



m\x + cr(m — mi), x > a, 
m x, \x\ < <J, 

mix — <r(m — mi), x < —a. 



(11) 
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For our choice of a = 15.6 and (3 — 25.58, mo = —8/7, mi = —5/7, the Chua system (jlOllip displays chaos (the 
emergence of a chaotic 'double scroll' attractor has been observed both in numerical simulations and in experiments 
17|). In the case of (jll[) . the function 4> a (x) has two discontinuity points, i.e., at ±<r. In the general case in which v 
discontinuity points are present, v independent gradient descent relations can be derived for each of the points and 
simultaneously integrated (together with other v corresponding auxiliary equations, analogous to Eq. (|5|)) in order 
to identify them all. Yet, in the case of the Chua system (TlT)|) , the problem can be simply formulated in only one 
unknown a (i.e., we rely on the fact that the two discontinuity points are symmetrical with respect to zero). 
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FIG. 2: (a) and (b) show the time evolutions of xi(t) (in black, thin line) and yi(t) (in grey, thick line) respectively at the 
beginning (t 6 [0,30]) and at the end (t € [470,500]) of the run (in (b) the two curves are superposed). Plot (c) shows 
|a:i(i) — J/i (t) | versus t g [0, 500]. As can be seen, synchronization is achieved for t approximately greater than 10 2 . (d) shows 
the time evolution of cr'(t) (in grey), which converges to the true value of a = 1 (the cr-th ordinate line is plotted as a black 
thin solid line in (d)). The true system obeys Eq. (|10|l . (a, /3, mo, mi) — (15.6, 25.58, —8/7, —5/7), with a — 1 and the model 
system obeys Eq. (|12jl . where the estimate a' evolves according to Eqs. (|14l llfcif) . /i(x) = xi, T = [7, 0, 0] T , 7 = 15, f = 1. The 
initial conditions for the true and the model systems are randomly chosen points from the Chua double scroll attractor, a' is 
evolved from an initial value that is far away from the true value of cr, i.e., er'(0) = 2.5, qi(0) = 0. 



We assume to model the true system by y(t) = [yi(t),y2(t),ts(t)] T , y(t) = F a >(y(t)), where a' is an estimate of the 
unknown true parameter a. We design an adaptive strategy to dynamically adjust a' to match the unknown value of 
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a. To this aim, we perform a one way diffusive coupling from the true system to the model, as follows, 

y(t)=^(y(t))+r[fc(x)-/i(y)], (12) 

where h(x) is in general an m < n vector of m observable scalar quantities that are assumed to be known functions 
of the system state x(£). r is an n x m constant coupling matrix. In what follows, we assume for simplicity that h is 
a scalar function (m = 1), h(x) = x\, T — [7, 0, 0] T and we proceed under the assumption that the value of 7 is such 
that when a' = a, the synchronized solution, y = x, is stable. We introduce the following potential, 

tt = [fc(x)-A(y)] 2 . (13) 

Again, W > by definition, and \& = if h(y) = h(x), that is, when the true system and the model system are 
synchronized. Thus, we seek to minimize the potential by adaptively evolving a' according to the following gradient 
descent relation, 

o' = -C0 = 2 C[Mx) - h(y)]Dh(y)^ = 2C[/i(x) - /i(y)]£>/i(y)q, (14) 

£ > 0. Note that for our choice of h, Dh = [1, 0, 0], yielding Dhq = qi, where q\ is the first component of the vector 
q. Therefore, we seek a differential equation that describes how qi = dyi/da' evolves in time. We note that (jlip can 
be rewritten as, 

<Pa(yi) = m yi + (m 1 - m )[H(-a' - yi)(yi + a') + H(y% - <j')[yi - a')]. (15) 

Then, from Eqs. (|12ll5p . we obtain the following differential equation for qi, 

qi(t) = c(t) qi (t) + d(t), where, (16a) 
c(t) = -7 - a{l + m Q + (mi - m )[H(-a' - Vl ) + H( yi - a')}}, (16b) 
d(t) = -a(ml ~ mO)[H(-a' - y t ) - H(yi - a')}. (16c) 



Our adaptive strategy is then fully described by Eqs. (Tj"2l \T%[ [To]) . In order to test the strategy, we run a numerical 
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experiment (shown in Fig. 2), in which we integrate Eqs. (| 101 121 [T4l [T6|) . We initialize both the true and the model 
systems from uniformly distributed random initial conditions on the Chua chaotic attractor, and a' is evolved from 
an initial value that is far away from the true value of a — 1, i.e., cr'(0) = 2.5, <?i(0) = 0. As can be seen, our adaptive 
strategy is successful in identifying a and in achieving synchronization between the true and the model systems. 

Hereafter, we consider a completely different identification problem, in which the dynamics of the system to be 
identified is described by a set of delay differential equations, 



where x = (xi, x 2 , x n ) T , $ S R n — > R n , r is a time delay. 

We try to model the dynamics of the true system by y(t) — 3>(y(i),y(i — r')), where t' is an estimate of the 
unknown true coefficient r. Our goal is to evolve r' in time to match the unknown true r-value and in so doing, to 
achieve synchronization between the model and the true systems. Therefore, we perform a one way diffusive coupling 
from the true system to the model, as follows, 



where h(x) and T are the same as defined before. In what follows, we assume for simplicity that h is a scalar function, 
i.e., m = 1. 

We now introduce the potential 'J, defined in (|13[) (similar considerations apply as in the previous cases) and we 
propose to minimize \P by making r' converge onto the true value t, through the following gradient descent relation, 



(3 > 0. We are now interested in how the n-vector r = dy' jdr' evolves in time. We note that the following equation 
describes the evolution of r(<), 




(17) 



y(t) = $(y(t), y(t - r'(t))) + T[h( X (t)) - h(y(t)% 



(18) 



P— = 20[h(y(t)) - h^mDhjfj = 2/3[h(y'(t)) - h(y(t))]Dhr, 



(19) 



*(*) - [ 



d<t> 



TDh]r(t) 



y'(t-r). 



(20) 



dy'{t) 



dy'(t - t) 



Therefore our adaptive strategy is fully described by the set of equations (|18I19I20[) . 
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We now present a numerical experiment testing the above strategy for a case in which the unknown true system 
is described by the Mackey- Glass equation, for which n = 1, y(t) = y(t), <p(y(t),y(t — r)) = —by(t) + ay(t — r)/(l + 
y(t — r) 10 ), T coincides with the scalar coupling 7. We choose a — 0.2, b — 0.1, r = 23, h(y) = y. For these values of 
the parameters, the dynamics of the Mackey-Glass system is chaotic. 

Both the true system and the model system are initialized from uniformly distributed random initial conditions 
on the chaotic attractor. We have preliminarily found that for r' = r = 23, the two systems synchronize for 
7 > 7* = 0.07. Thus, for our numerical experiment, we choose 7 = 0.1. Moreover, we take = 1. The results of our 
numerical experiment are shown in Fig. 3. The experiment consists of two parts. For t < 10 4 , r is kept constant and 
equal to 23, while r' is initialized from a value that is far away from the true value of r = 23, i.e., r'(0) = 15. As 
can be seen, after a transient, t >0.25 x 10 4 , r' converges to the true value of r = 23 and the synchronization error 
\x — y\ decays to zero. For 10 4 < t < 2 x 10 4 , r becomes a function of time, i.e., r(t) = 23 + 3 sin(27rl0~ 4 t), while the 
adaptive strategy (|18ll9l20p is kept running. As can be seen, for t > 10 4 , r'(t) tracks quite well the time evolution of 
r(t) and approximate synchronization between the model and the true system is attained. 

Synchronization has been shown to be a convenient tool to identifying the dynamics of unknown systems. Different 
techniques have been proposed, see e.g., . Yet, the problems of identification of delays and discontinuity points 

that commonly characterize the dynamics of real systems have not received adequate attention in the literature. In 
this paper, we have proposed a general methodology that, based on a simple gradient descent technique (MIT rule), 
can be successfully employed to address both these problems. We have also shown the usefulness of our strategy in 
the case that the unknown parameters of the true system to be estimated slowly drift in time (here, by 'slowly' we 
mean that they evolve on a time scale which is much longer than that on which a typical chaotic oscillation occurs). 

We successfully apply our methodology to solve problems as different as the identification of delays and the identi- 
fication of discontinuity points of unknown dynamical systems. This suggests that our approach can provide a simple 
and flexible paradigm for the resolution of a variety of diverse identification problems. Furthermore, by making use 
of the properties of synchronization of chaos, the adaptive strategy can be used to extract several unknowns from 
knowledge of the chaotic state-evolution of the system of interest. 
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